clc;
clear;

basis_type = "Hermite";
derivative_type = "x";
process_showbasis(basis_type, derivative_type);





%% process_showbasis
function process_showbasis(basis_type, derivative_type, varargin)
[x, y] = meshgrid(0:0.0125:1, 0:0.0125:1);
z = basis_function_ref(x(:), y(:), basis_type, derivative_type);
if isempty(varargin)
    f_begin = 1;
    f_end = size(z,2);
else
    f_begin = varargin{1};
    f_end = varargin{1};
    if f_end > size(z,2)
        f_end = size(z,2);
    end
end
triangle_region = @(x,y) (x >= 0) & (x <= 1) & (y >= 0) & (y <= 1) & (y <= -x+1);
for j = f_begin:f_end
    Z = reshape(z(:,j), size(x,1), size(x,2));
    Z(~triangle_region(x, y)) = NaN;
    surf(x, y, Z);
    hold on
end
end

%% basis_function_ref
function N = basis_function_ref(x, y, basis_type, derivative_type)
if strcmp(basis_type, "Non-com")
    N = zeros(size(x,1), 3);
elseif strcmp(basis_type, "Hermite")
    N = zeros(size(x,1), 10);
else
    error("Invalid basis type.");
end
if strcmp(basis_type, "Non-com")
    if strcmp(derivative_type, "x")
        N(:,1) = 2*x+2*y-1;
        N(:,2) = 1-2*x;
        N(:,3) = 1-2*y;
    elseif strcmp(derivative_type, "dx")
        N(:,1) = 2;
        N(:,2) = -2;
        N(:,3) = 0;
    elseif strcmp(derivative_type, "dy")
        N(:,1) = 2;
        N(:,2) = 0;
        N(:,3) = -2;
    end
elseif strcmp(basis_type, "Hermite")
    lambda = 1-x-y;
    if strcmp(derivative_type, "x")
        a = x.*y.*lambda;
        N(:,1) = lambda.^2.*(3-2*lambda)-7*a;
        N(:,2) = x.*lambda.^2-a;
        N(:,3) = y.*lambda.^2-a;
        N(:,4) = x.^2.*(3-2*x)-7*a;
        N(:,5) = x.^2.*(x-1)+2*a;
        N(:,6) = x.^2.*y-a;
        N(:,7) = y.^2.*(3-2*y)-7*a;
        N(:,8) = y.^2.*x-a;
        N(:,9) = y.^2.*(y-1)+2*a;
        N(:,10) = 27*a;
    elseif strcmp(derivative_type, "dx")
        b = y.*(lambda-x);
        N(:,1) = 6*lambda.*(lambda-1)-7*b;
        N(:,2) = lambda.*(lambda-2*x)-b;
        N(:,3) = -2*lambda.*y-b;
        N(:,4) = 6*x.*(1-x)-7*b;
        N(:,5) = x.*(3*x-2)+2*b;
        N(:,6) = 2*x.*y-b;
        N(:,7) = -7*b;
        N(:,8) = y.^2-b;
        N(:,9) = 2*b;
        N(:,10) = 27*b;
    elseif strcmp(derivative_type, "dy")
        c = x.*(lambda-y);
        N(:,1) = 6*lambda.*(lambda-1)-7*c;
        N(:,2) = -2*x.*lambda-c;
        N(:,3) = lambda.*(lambda-2*y)-c;
        N(:,4) = -7*c;
        N(:,5) = 2*c;
        N(:,6) = x.^2-c;
        N(:,7) = 6*y.*(1-y)-7*c;
        N(:,8) = 2*x.*y-c;
        N(:,9) = y.*(3*y-2)+2*c;
        N(:,10) = 27*c;
    end
end
end